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ABSTRACT 


A  unified  approach  to  the  dynamic  nonlinear  estimation  problem  is 
presented.  General  algorithms  for  filtering,  smoothing,  and  prediction  for 
both  continuous  and  discrete  nonlinear  systems  are  derived.  The  continuous 
problem  is  treated  first  by  a  direct  approach  which  yields  an  exact  solution. 

A  similar  approach  is  then  applied  to  the  more  complex  discrete  estimation 
problem.  The  discussion  provided  points  out  some  subtle,  but  important 
features  of  the  main  results  which  should  be  useful  when  applying  the  results 
to  practical  problems.  To  provide  insight,  motivation  and  background  for  the 
technique  used  in  obtaining  the  solution  of  the  general  nonlinear  problem,  the 
special  case  of  linear  estimation  is  discussed  in  the  first  two  appendices.  In 
a  third  appendix  it  is  shown  that  in  the  continuous  measurement  limit,  the 
derived  discrete  and  continuous  nonlinear  filtering  algorithms  are  identical. 
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Discrete  Nonlinear  Estimation 


I.  INTRODUCTION 


For  nearly  two  centuries  the  concept  of  least  squares  has  served 
scientists  well  in  the  solution  of  estimation  problems.  Although  the  technique 
is  basically  linear,  it  has  been  used  effectively  in  nonlinear  problems  by  suc¬ 
cessive  application  to  linearized  approximations  (differential  corrections). 

In  recent  years,  however,  the  requirements  for  increased  accuracy  and 
efficiency,  particularly  in  the  aerospace  industry,  have  created  considerable 
interest  in  the  development  of  a  purely  nonlinear  technique.  With  few  excep¬ 
tions,  the  least  squares  concept  has  been  ignored  and  probability  theory  has 
been  adopted  as  the  basis  for  new  developments.  Although  this  recent  approach 
has  yielded  some  interesting  mathematics,  the  results  reported  to  date  do  not 
give  a  practical  complete  approach  to  nonlinear  estimation  (3ee,  for  example, 
comments  by  R.  S.  Bucy  [*]).  In  any  case,  with  the  use  of  probability  theory 
there  is  the  inherent  assumption  that  a  system's  dynamics,  measurements, 
and  relevant  statistics  can  be  modeled  perfectly.  Clearly  this  is  not  generally 
the  case  in  practical  engineering  situations.  This  report  will  show  that  in 
abandoning  the  least  squares  approach,  many  researchers  have  overlooked  a 
very  practical  and  concise  nonlinear  estimation  technique  which  is  not  limited 
by  modeling  assumptions. 

An  important  breakthrough  in  extending  least  squares  directly  to  non¬ 
linear  problems  was  made  by  Detchmendy  and  Sridhar  in  1965  [5]  ■  Their 
results,  however,  are  limited  to  filtering  for  continuous  nonlinear  .systems, 
and  even  then  their  derivation  is  indirect  and  unnecessarily  involved. 

In  the  first  part  of  this  report  the  continuous  filtering  problem  will  be 
solved  from  a  new  viewpoint.  A  similar  approach  will  then  be  used  to  derive 
an  appropriate  filter  for  discrete  nonlinear  systems.  Algorithms  for  smooth¬ 
ing  in  both  continuous  and  discrete  nonlinear  systems  will  also  be  presented. 
The  discussions  provided  point  out  some  subtle,  but  important,  features  of 
the  overall  technique.  They  should  be  helpful  in  application  of  the  main  results 
to  practical  problems.  Three  appendices  are  also  included.  The  first  two 
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furnish  the  necessary  background  and  motivation  for  the  technique  used  in 
solving  the  nonlinear  estimation  problem.  The  last  one  shows  that  in  the 
continuous  measurement  limit  the  derived  discrete  and  continuous  nonlinear 
filtering  formulas  are  identical. 

BASIC  NOTATIONAL  CONVENTIONS 


a)  Matrices  are  upper  case  letters. 

b)  Vectors  are  lower  case  letters. 

c)  Exceptions  to  these  rules  are  i,  k,  n,  p,  q  which  are 
used  as  integers;  t,  T  which  denote  time;  and  J  which 
represents  a  scalar. 

d)  Asterisks  (*)  denote  matrix  transposition. 

e)  Vectors  are  assumed  to  be  columns  unless  otherwise 
denoted  by  an  asterisk. 

f)  Parentheses  (  )  are  used  exclusively  for  showing 
functional  dependence,  while  square  brackets  [  J 
are  used  for  grouping  terms  in  equations. 

g)  The  operator  d/dx  applied  to  a  scalar  (vector,  matrix) 
gives  the  appropriate  gradient  (Jacobian,  Hessian) 
evaluated  at  x. 


II.  CONTINUOUS  STATE  ESTIMATION 


A.  DERIVATION  -  CONTINUOUS  SYSTEMS 
Consider  the  system 


x  =  f  (x,  t)  +r(t)w(t) 

(1) 

z  =  h(x,  t)  +  v(t) 

(2) 

where  x  is  a  p-dimensional  state  vector,  z  is  a  q-dimensional  vector  of 
measurements,  f(x,t)  and  h(x,t)  are  appropriately  dimensioned  vector  functions 
of  x  and  t,  and  r(t)w(t)  and  v(t)  are  vectors  representing  unknown  disturbances, 
noise,  and  general  modeling  errors,  with  T(t)  being  a  known  px  k  matrix 
function  of  time.  It  is  to  be  emphasized  that  w(t)  and  v(t)  are  not  necessarily 
random  functions.  They  will,  however,  be  assumed  to  have  a  zero  time  average. 
Note  that  a  non-zero  time  average  can  be  estimated  as  part  of  the  state  vector. 

The  generalized  least  squares  estimation  problem  can  be  stated  suc¬ 
cinctly  as  follows.  Find  the  y  and  u  which  minimize  the  scalar 


J  =  2  "  x0]  P0  [yV  "  X0] 

X 

z  -  h(y ,  t) J  R-1(t)[z  -  h( y,  t) J  +u  Q_1(t)uj  dt 


(3) 


subject  to  the  constraint 


y  =  f(y,t)  +  r(t)u(t)  (4) 

where  Pq,  R,  and  Q  are  arbitrary  symmetric  positive  definite  weighting 
matrices  and  Xq  is  an  a  priori  guess  of  the  initial  state  of  the  system.  As 
stated,  the  problem  is  strictly  one  of  smoothing;  i.e.,  the  solution  at  any 


1/  WtMMFmmmmnmtmmrn 


instant  t  will  be  the  least  squares  estimate  of  the  state  vector  based  on  the 
data  gathered  in  the  entire  interval  t^  £  t  £T,  The  filtering  aspects  of  the 
problem  will  become  apparent  later. 

Note  from  the  form  of  (3)  and  (4)  that  the  problem  can  be  considered 
purely  as  one  in  control  theory;  i.e.,  it  can  be  considered  one  of  finding  a 
control  u(t)  to  force  the  output  of  the  given  system  to  track  an  arbitrary  signal 
z(t).  This  viewpoint  is  particularly  useful  for  indicating  the  possible  applica¬ 
tion  of  the  function  u(t)  to  the  modeling  problem.  Namely,  since  u(t)  is  the 
"best  control"  for  driving  a  system  with  dynamics  f(y,t)  to  track  z(t),  it 
can  be  considered  an  estimate  of  the  difference  between  f(y,t)  and  the  true 
dynamics  of  the  system  generating  z(t). 

The  given  problem  of  (3)  and  (4)  clearly  calls  for  application  of  the 
calculus  of  variations.  By  a  well-known  procedure  [1,  3,  6]  it  can  be  shown 
that  the  solution  must  satisfy  the  two-point  boundary  value  problem  (TPBV) 

y  =  f(y.t)  -  TQr\  (5) 

^  x  +  |^  R"lU  -  Wy.t)]  (6) 

x<v =  ■  {MT)  =  0)  (7) 

with 

u(t)  =  -  Qr*\(t)  (8) 

where  X  is  a  p- vector  of  unknown  functions  (adjoint  variables).  In  general, 
a  direct  solution  of  an  arbitrary  TPBV  problem  is  quite  unlikely.  For  this 
case,  however,  a  special  approach  leads  to  a  very  practical  method  of  gen¬ 
erating  solutions.  The  procedure  used  is  motivated  by  the  solution  of  the 
linear  version  of  this  problem  for  which  details  are  given  in  Appendix  A. 
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Assuming  continuity  of  the  functions  involved,  these  equations  are  then 
necessary  and  sufficient  to  ensure  the  validity  of  (9)  in  a  neighborhood  of  \  =  0. 
Consequently,  by  (7)  they  must  hold  near  t  =  T,  the  end  of  the  data  span.  Fur¬ 
thermore,  the  term  "filter"  is  correct  for  (13)  and  (14)  since  by  (9),  X  =  0 
gives  y  =  x,  the  best  estimate  of  the  state  vector  at  T  based  on  data  gathered 
in  the  interval  t^S  t  i  T. 

To  obtain  initial  conditions  for  these  equations,  observe  that  T  can  be 
considered  arbitrary,  so  that  if  T  —  tg,  it  must  follow  that 

x(tQ)  =  xQ  (15) 

p(tQ)  =  P0  (16) 

Note  also  from  (9)  and  (7)  that  x(T)  is  the  missing  boundary  condition 
for  (5)  so  that  a  backward  solution  of  (5)- (6)  using 

\(T)  =  0  and  y(T)  =  x(T)  (17) 

immediately  gives  the  complete  smoothing  solution  in  the  interval  tg  £  t  £  T. 

These  results  are  summarized  in  Table  I  using  the  standard  notation 
x(t/T)  for  the  state  estimate  at  time  t  based  on  data  up  through  time  T.  For 
completeness  the  solution  to  the  prediction  problem  is  presented. 

B.  DISCUSSION  -  CONTINUOUS  SYSTEMS 

The  basic  formulas  for  filtering  (18)  and  (19)  have  been  derived  pre¬ 
viously  by  Detchmendy  and  Sridhar  [5] .  Their  derivation,  however,  is 
incomplete  because  it  does  not  yield  the  all-important  initial  conditions.  In 
addition,  their  approach  is  unnecessarily  complicated  by  the  introduction  of 
a  partial  differential  equation. 

A  potential  application  of  least  squares  estimation,  which  apparently 
has  been  overlooked,  is  its  possible  use  in  the  solution  of  the  modeling 
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Table  1.  Continuous  Nonlinear  Estimation 


Filtering:  t  =  T 


dx 

cTT 


f(x,t)  +  P  ^-R’^z  -  h(x,t)] 
dx 


(18) 


dP 

ar 


—  P  +  p  +  ror*  +  P 

a£  i* 


R-1[ z  -h(£,t)])p  (19) 

dx  \  dx  / 


*V  =  xQ  and  P(tQ)  =  PQ 


(20) 


Smoothing:  t  <  T 


d\ 

dt 


^^(t/T)  =  f(A(t/T),t)  -  ror\ 


ai*  *  .  dh*  „-i 


dx(t/T) 


\  + 


dx(t/T) 


R’  [z  -  h(x(t/T),  t)] 


(21) 

(22) 


x(T/T)  =  x(T)  and  X(T)  =  0 


(23) 


Prediction:  t  >  T 


^S(«/T) 


=  f(x(t/T),  t) 


x(T/T) 


=  x(T) 


(24) 

(25) 


..  *4|> 


if*  VM 


mm**  «»ip«wm 


problem.  At  present,  errors  in  modeling  system  dynamics  are  often 
the  most  prominent  cn.use  of  inaccuracy  in  practical  estimation  problems. 

In  the  least  squares  formulation  of  the  problem,  (l)-(4),  the  general  system 
modeling  error  is  included  in  the  term  w.  The  solution  in  Table  1  offers 
both  a  filtered  and  a  smoothed  estimate  of  w.  From  (18)  and  (8)  these 
are 

Iw(t)  =  P-^  R"l[z  -  h(£,t)]  (26) 

d£ 


and 


r&t/T)  =  -  ror*\(t)  (t  <  T)  (27) 

respectively.  Note  that  while  the  solution  to  the  prediction  problem  (24) 
follows  from  (5)-(6)  with  z(t)  =  h(y,t)  =  0,  it  may  not  be  the  most  practical 
one  for  modeling  problems.  That  is, if  rwis  significant,  it  may  be  advan¬ 
tageous  to  extrapolate  it  in  some  consistent  manner  for  use  in  state  vector 
prediction.  One  approach  would  be  to  find  an  appropriate  dynamic  model  for 
this  error,  then  augment  the  given  state  vector  accordingly,  and  apply  the 
formal  prediction  algorithm  (24). 


III.  DISCRETE  STATE  ESTIMATION 


A.  DERIVATION -DISCRETE  SYSTEMS 

Although  the  continuous  estimation  formulas  of  Table  I  can  be  useful  in 
certain  situations,  there  are  many  practical  estimation  problems  in  which 
system  measurements  are  not  recorded  continuously,  but  rather  are  known 
only  at  specific  instants  of  time.  This,  coupled  with  the  fact  that  it  is  generally 
desirable  to  perform  analyses  digitally,  is  good  reason  to  seek  a  discrete 
analogue  of  the  least  squares  formulas  of  Table  I.  Several  attempts  at 
deriving  such  formulas  have  been  made  in  the  past  [11,  12]  with  no  overall 
success.  This  is  not  surprising  since  the  extension  is  not  trivial  as  the  dis¬ 
cussion  below  will  show.  Although  the  discrete  problem  is  similar  to  the 
previous  one,  there  are  sufficient  differences  to  warrant  a  completely 
independent  treatment.  Since  the  results  to  be  obtained  are  new,  a  par¬ 
ticularly  detailed  discussion  is  provided. 

For  the  discrete  problem  consider  the  system 

xi  =  gi-l<xi-l*  +  ri-lwi-l  (2  s  i  £  n)  (28) 

z.  =  h.lx^  +  vi  (l  <  i  <  n)  (29) 

where  x.  is  a  p-dimensional  state  vector,  z.  is  a  q-dimensional  vector  of 
i  i 

measurements,  g._j(x._j)  and  lv(x^)  are  appropriately  dimensioned  vector 

functions,  and  I\  jw._j  and  v^  are  vectors  representing  unknown  disturbances, 

noise,  and  general  modeling  errors,  with  T.  ,  being  a  known  p  X  k  matrix 

function  of  time.  It  is  to  be  emphasized  thatw^  ^  and  are  not  necessarily 

random  functions.  They  will,  however,  be  assumed  to  have  a  zero  time  average. 

Note  that  a  non-zero  time  average  can  be  estimated  as  part  of  the  state  vector. 

In  this  context  the  generalized  least  squares  problem  is  that  of  finding 

the  sequences  y,iy,,  .  .  .  ,  y  and  u4,u, ,  .  .  .  ,  u  .  which  minimize  the 
i  c  n  14  n-i 

scalar 
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J  =7tyl  ■  *01#P0ltyl  “  *03  +  Itzl  "  hl(yl)1*Ri 1[ Z1  •  hl(yi)] 
n 

+  ?  E  [k  -  WlX ‘K  •  hit',i)1  +  ui-iQi-iui-t] 


subject  to  the  constraints 


yi  =  *i-i(yi-il  +  ri-.Vi 


(2  £  i  <n) 


where  P^,  R.  and  ^  are  arbitrary  symmetric  positive  definite  weighting 
matrices,  and  Xq  is  an  a  priori  guess  of  the  initial  state  of  the  system  at 
time  tj .  Observe  that  (28)  is  meaningless  for  n  =  1,  and  the  basic  problem 
reduces  to  finding  the  minimum  of 

J1  =Ityl  -*01*P01(yl  ’*0]  +Itzi  -  h1(y1)]*R[1[zl  -  h1(y1>]  (32) 


The  multistage  problem  of  (30)  and  (31)  clearly  falls  within  the  scope 
of  the  well-known  Lagrange  multiplier  technique  for  solving  constrained 
extrema  problems.  By  introducing  two  dummy  constraints,  y^  =8q(Vq)  +  uq 
and  yn+j  =  g^(yn>  +  un»  a  general  formulation  can  be  found  which  conveniently 
includes  the  single  stage  problem  (32).  In  particular,  write 

J'  =7^1  •i0)*Po‘lri  ‘40J  +  Cllyn+1  -  8„<yn»  -  rnun]  +iuOQo‘uO 


n 

+  -  h.(y.)]*R:1[I.  -  h.(y.)l  +KQi‘ui 

+  xi*[yi‘gi-t(yi-J)-ri-lui-l]] 
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where  \^,X2»  ....  are  p- vectors  of  Lagrange  multipliers.  Now 

both  J  and  Jj  will  be  minimized  when  8J'/8y^  =  0  for  i  =  0, 1 , 2,  .  .  .  ,  n  +  1 ; 
and  8J'/8u.  =  0  for  i  =  0,  1 , 2 ,  .  .  .  ,  n,  or,  more  explicitly,  when 


dgi  A 

X1  =  "aTj-  X2  +  9y”  R1  P  "  hl(yl^  "  P0  t  yl  '  x0] 


M 

dgi 


dh. 


X.  =  ~-X..  +  -jP-  r:1  [z.  -  h  ( y  ) ] 

l  8y^  it  l  oy^  l  i  li 


y.  =  g.  ,(y.  . )  +  r.  .Q.  .  r.  .  X. 

7i  °i-l  7i-l  i-I  i-I  i-1  i 


xi'W> 


(2  <  i  <  n) 


(2  <  i  <  n) 


(34) 


(35) 


(36) 


(37) 


with 


u0  =0 


and 


Ui  =  QiFi  \+l 


(1  <  i  S  n) 


(38) 


Note  that  (34),  (37)  and  (38)  show  Jj  is  minimized  for  n  =  1. 

Observe  that  the  problem  described  by  (34)  through  (37)  is  analogous 
to  the  TPBV  problem  encountered  earlier  in  the  continuous  estimation  prob¬ 
lem.  Recall  in  that  case  that  the  method  of  solution  was  motivated  by  an 
approach  to  the  corresponding  linear  problem.  A  similar  technique  can  be 
used  here.  A  solution  to  the  linear  discrete  estimation  problem  is  presented 
in  Appendix  B.  That  derivation  strongly  suggests  that  a  substitution  of  the 
form 


9gi 

y.  =  p.  P  x., .  +  x. 

1  IgA  1+1  1 

1 


(39) 
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d 


•Jfp  > 


„  88i.K  a 

ri-i  =pi-iX-Ki  +  xi-i 
^i-l 


(40) 


where  and  j  are  arbitrary  and  x^  and  x^  ^  are  filtered  state  estimates, 
be  applied  to  the  TPBV  problem  (34)-(37). 

Proceeding  along  these  lines,  then,  combine  (36)  and  (39)  to  obtain 


*  dg- 

T.  ,Q.  .T.  A.  =  P.  ~-X..  .  +  x.  -  g.  . (y.  .) 

l-i  1-1  1-1  1  1  qA  i+1  i  ®i-l  7i-l 

i 


(41) 


Use  (35)  to  eliminate  X^,  and  expand  gj_j(yj_j)  about  x._j.  Then 


*g{ 


dh. 


ri- iQi. l ri- 1 "5^"  xi+  i  +  ri.A.,rr.<  Ri  1[zi  ■  hi(yi)J  = 

„  ^®i  %  .A  /A  i  ®®i-l  r  A  1 

P.  — —  X.  ,  |  +  X.  -  g.  .  X.  .  )  -  ■  ■  ly.  ,  -  X.  .  I 

1  „A  1+1  1  61-1  l-r  -A  lyl-l  l-lJ 

1  1- 1 


(42) 


At  this  point  the  index  i  is  confined  to  be  near  n,  i.e. ,  i  =  n  -  1,  n, 

where  n  is  the  last  stage  of  the  estimation  problem.  Other  values  of  i  are  of 

no  importance  in  the  discussion  below.  With  this  in  mind,  it  is  assumed  that 

y  ,  -  x  ,  and  X  are  "small"  in  the  sense  that  quadratic  and  higher  order 
'n-4  n-c  n 

terms  in  these  quantities  can  be  ignored.  The  implications  of  this  hypothesis 
will  be  considered  later.  Note  from  (40)  that  assuming  to  be  small  guar¬ 
antees  that  y  ,  -  x  4  will  also  be  small.  Thus,  for  the  i  of  interest, 
n*  i  n*  •  ^ 

terms  nonlinear  in  y.  ^  -  x.  j  in  (42)  can  be  omitted. 


Continuing  with  the  development,  then,  apply  (40)  and  (35)  to  (42)  to  get 


...  dh. 


ri-iQi-iri-.  ir-xi+i  +  ri.iQi.iri*-i  WJ  * 


^  *i  »  A  a  v 

P.  — —  \..  .  +  X.  -  g.  .  (x.  .) 
1  -A  1+1  1  61-1'  1-1 

dx. 
i 


.  ("  3g-:  ,  rag:,!  dh.  .  1 

p-  i  a-1"  x  +  a^-  R"'[2-  -  h.fy.)J 

.  ,  ‘-1 3$.  ,  L0yi  1+1  9n  1  1  1  1 J 

l-i  L  i-i  j 


Now,  bearing  in  mind  the  fact  that  yi  is  a  function  of  \.+  1  by  (39),  and  that 
V+j  is  assumed  to  be  small,  expand  both  sides  of  (43)  in  powers  of  to 
obtain 


dh.  .  a  / dh .  •  ,  \  dg. 

-  h^‘,i  +  ifcr  *»  f 


i  i 


*  a8i  a8i  A  A 

+  r.  ,Q.  ,r.  .  =  p.  .  +  x.  -  g.  ,(x.  .) 

1-1  1-1  1-1  aA  l+i  1  A  1+1  1  sl-l  1-1 

dx.  dx. 

1  1 


■  ,  dg.  ,  f  dg.  dh.  , 

i •  •  n  1"  i  I  1  \  i  i  «  “i  r  i  /A  v  i 

-  Pj  1  “T — 7”  A...,  +  — T—  R.  Iz.  -h.  x.  J 

i_.  L  A-.K  4  11  11 

8«.W  11  1  1  /  >a£.  1+1 


For  convenience,  the  matrix  is  defined  as 


A  *  °»i-l  °Si_l 

pi/i-i  ■  ri.iQi.iri.i ♦  pr“pi-i  7^ 

3xi-i 
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and  the  vector  as 


lify*  T'l'l'i-bA11 

®*i 


(46) 


Then  (44)  can  be  expressed  more  concisely  as 


8gi  A  d\  88;  ag; 

p. ,  —■*  X..  .  +  P.,.  ^.(X.)  fP.,.  .  “^P.  -JT-X..  4  =  P.  “ rr“X  . ,  . 
l'l“l  1+1  1/1-1  1  1  1/1-1  aA  1  i+l  1  QA  l+l 


axi 


ex.  ax. 


axi 


i  A  #A  i 

+  x.  -  g.  ,  (x.  .) 
1  ®i-i  i-r 


(47) 


Let  i  =  n.  Then  \n+J  =  0  by  (37),  and  (47)  becomes 


x  =  g  t  (x  * )  +  P  ,  |  r|  (x  ) 

n  °n-l  n-1  n/n-i'n  n 


(48) 


Now,  assume  that  new  data  is  acquired,  i.e.,  the  last  stage  becomes 
n  +  1.  Then,  Xn+£  =  0,  but  \r+j  ^  0,  so  that  equating  coefficients  of  X^+1  in 
(47)  for  i  =  n  gives 


0Tln 

Pn  =  Pn/n-l  +  Pn/n-l  ^  Pn 

ii 


(49) 


Thus  (45)  and  (46)  with  i  =  n,  and  (48)  and  (49)  are  a  set  of  recursion  relations 

A 

Note  that  (48)  '  limited  from  a  practical  viewpoint  because  it 


for  P  and  x  . 
n  n . 


is  implicit  in  xr.  Fortunately,  the  relationship  can  be  solved  directly  by 
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means  of  some  interesting  manipulations.  Consider  the  quantity  -  xn  where 
xn  ^  gn-l^*n-l^*  From  (36)»  O7)  311111  (40) 


-  x„  *  +  rn-lQn-irn-lXn  ’ 


Expanding  the  first  term  on  the  right  about  and  applying  (40)  and  (45) 

yields 


A  —  . 

X  -  X  =  P  /  .  X 

n  n  n/n-1  n 


where  X  has  already  been  assumed  to  be  small.  Thus  an  expansion  of  (48) 

^  A  — 

in  powers  of  x  -  x  can  be  truncated  after  the  linear  term  with  no  loss  of 
r  n  n 

generality.  In  particular,  (48)  becomes 


&  =  x  +  P  >  . 

n  n  n/n-1 


f  .  .  %  r  A  -  ,1 

1n<xn)+  lxn-xJ 

L  j 


or,  finally 


A  _  f  dq  "1-1 

X  =X  +  I  -  P  ,  .  - —  P  ,  .  T|  (x  ) 

n  n  n/n-1  ^  n/n-l'n  n 

■  M 


where  the  indicated  matrix  inverse  is  assumed  to  exist.  A  practical  inter¬ 
pretation  of  this  hypothesis  will  be  given  later. 

To  be  consistent  with  the  above  result  (49)  must  also  be  expanded  in 

powers  of  x  -  x  .  At  first  it  appears  that  third  derivatives  of  h  (x  )  will  be 
*  n  n  n  n 

involved  in  the  coefficients  of  the  linear  term.  Recall,  however,  that  (49) 
resulted  from  equating  coefficients  of  X^j  in  (47).  Since  that  quantity  itself 


is  assumed  to  be  small,  inclusion  of  the  first  power  of  x  -  x  would  then 

n  n 

effectively  result  in  a'feecond  order"term.  Thus,  (49)  immediately  reduces  to 


\>=  [1-pn/n-l  'Pn/n-l 

L  °n  J 


Note  that  (53)  can  now  be  written  as 


x  =x  +  Pr](x) 
n  n  n  'n  n 


Thus  (54),  (55)  and 


*  8*n-l  8gn-l 

Pn/n-l  2  rn-lQn-irn-l  +  Pn-1 

°n- 1  °n- 1 


form  the  desired  set  of  recursion  relations  for  filtering  nonlinear  systems. 

To  establish  initial  conditions  for  the  sequence,  let  n  =  2  and  consider 
the  original TPBV  problem  (34)-(37).  From  (34)  with  =0,  it  is  clear  that  y j  is 
some  function  of  \2<  Expanding  y^  in  powers  of  gives  y^  =  +  [dy^/d^]^ 

+  .  .  .  ,  where  x^  must  be  the  solution  to  the  single  stage  estimation  problem. 
Since  by  hypothesis  is  small,  this  result  shows  that  y^  -  x^  must  also  be 
small.  Applying  these  ideas  to  an  expansion  of  (34)  in  powers  of  yt  -  x^ 
yields 


8gT  8r)| 

0=3rx2  +  r,i(xi)+-7"[yi  -xi] 


P0  f  xl  ■  x0)  ■  pot  yi  ■  XJ 
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Since  by  definition  Xj  satisfies 

-  *0i  -° 

(57)  reduces  to 


(58) 


(59) 


Now  assume  that  -  x^  is  small  and  expand  Hjfxj)  and  dr^/Ox^  about  Xg. 
Thus  (58)  and  (59)  become 


(60) 


(61) 


(62) 


Thus  (60)  and  (62)  are  the  required  initial  conditions  for  the  filter.  Note  that 
the  assumption  that  -  Xg  is  small  is  totally  consistent  with  those  made 
earlier  in  developing  the  basic  recursion  relationships. 

As  in  the  continuous  case,  the  solution  to  the  smoothing  problem  follows 
directly  from  the  results  of  filtering.  That  is,  xr  and  ^-n+j  =  0  are  final  con¬ 
ditions  for  running  (35)  and  (36)  backward  from  i  =  n. 
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(63) 

(64) 

(65) 

(66) 

(67) 

(68) 

(69) 

(70) 

(71) 

(72) 
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For  the  convenience  of  the  reader,  the  results  for  the  discrete  estimation 
problem  are  summarized  in  Table  II.  The  notation  x^n  is  used  to  indicate  an 
estimate  of  the  state  vector  at  t.  based  on  data  through  tn-  For  completeness 
the  prediction  results  are  included.  It  is  noted  that  in  the  filter  is  itself  a 
prediction  which  is  updated  by  the  measurement  z^  to  give  x.. 

Observe  that  the  use  of  (64)  and  (65)  requires  inversion  of  a  p  X  p 
matrix  at  each  stage.  If  3(9h. /0xi)/3xi  or  z.  -  h.  (x^  is  small  (this  includes 
the  case  when  measurements  are  linear),  then 

/  8hV  .  \  9h"  .  9h. 

£(  £ 


and  the  famous  Matrix  Inversion  Lemma  [9]  can  be  applied  to  (64)  and  (65) 
to  yield 


A  _  ah'."  9h.  9h“ 

x.  =  x.  +  P.  .  — —  — —  P.  ,.  .  — -  +  R.  [z.  -  h.(x.)]  (74) 

1  1  i/1-1  a-  ar  l/i-1  a-  1  1  11 


9x.  9x. 

l  l 


9h?’  I  9h.  9h*  l'1  ah. 

P.  =  P.  -  P.  .  — -  — -  P.  4  — -  +  R.  — -P./. 


i  i/i-l  i/i-l  3x  9x  l/i-l  dx  i  9x  l/i-l 

ii  i  i 


These  formulas  require  inversion  of  q  x  q  matrices.  The  result  (75)  can,  in 
turn,  be  put  in  the  more  computationally  attractive  form 


9h.  3h.  3h* 

K.  =  P.  ,.  .  — -  — L  P.  ..  .  — -  +  R. 

1  l/l-‘  9x.  ax.  l/l-1  ax.  1 

1  1  1 


r  ahi  i  r  ahii* 

P.  =  I  -  K.  —  P. .  I  -  K.  —  +  K.R.K. 

1  L  1  ax.  J  l/l-!L  1  ax.  J  111 

u  1  w  1 


i 


The  consistency  of  the  approach  taken  here  with  that  used  earlier  for 
continuous  systems  is  verified  in  Appendix  C.  There  it  is  shown  that  the 
limiting  case  of  the  discrete  filters  of  Table  II  as  the  time  between  measure¬ 
ments  becomes  small  is  indeed  the  continuous  filter  of  Table  I. 

B.  DISCUSSION -DISCRETE  SYSTEMS 

1.  Fundamental  Interpretations 

In  this  section  the  assumptions  made  in  deriving  the  formulas  of  Table  II 
will  be  discussed  in  detail  along  with  other  items  of  general  practical  interest. 

First  consider  the  assumption  that  \  and  yn_2  -  xn_2  are  small.  From 
(35)  it  can  be  seen  that  X  will  be  small  if  £  3yn]R-n*[zn  “  hn(yn)]  *s  sma^‘ 
This  will  be  true  if  a  least  squares  fit  of  the  last  stage  data  alone  gives  results 
near  those  of  the  composite  n- stage  fit.  Thus,  small  X  (  is  essentially  a  require¬ 
ment  that  new  data  be  consistent  with  previously  processed  information.  Simi¬ 
larly,  the  assumption  of  small  yn2  -  xn_2  (i.e.,  closeness  of  filtering  and 
smoothing  near  n)  is  a  requirement  that  the  introduction  of  new  data  does  not 
significantly  change  the  results.  Note  by  comparison  to  the  derivation  of  the 
continuous  formulas  of  Table  I,  that  the  above  assumptions  effectively  simulate 
a  "continuity'1  condition  for  the  discrete  problem.  Basically  then, the  restrictions 
are  mathematical  in  nature  and  are  not  necessarily  limiting  from  an  engineer¬ 
ing  point  of  view.  Certainly,  if  difficulty  were  to  occur  in  a  practical  problem, 
it  would  be  expected  at  the  initial  stages  of  filtering.  An  approach  to  be  fol¬ 
lowed  in  such  situations  will  be  suggested  in  a  later  section. 

In  the  derivation  of  (63)  -  (66)  the  existence  of  the  indicated  matrix 
inverse  was  assumed.  This  can  be  justified  on  practical  grounds  as  follows. 

Let  i  =  1  and  note  that  Pj/q  is  positive  definite.  Then 


(78) 
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Now,  clearly  [  ahj/SXj  ]  Rj  3hj/3Xj  J  is  at  least  positive  semidefinite.  The 

matrix  sum  on  the  right  of  (78)  will  thus  be  positive  definite  unless  the  quantity 

Rj  *  [zj  -  hj(Xj)]  is  "large".  Such  a  situation,  however,  implies  that  the 

data  z.  is  especially  inconsistent  with  x.  .  It  could  then  be  argued  that  the 
1  - 1  1 

weighting  matrix  R^  be  decreased  to  limit  the  effect  of  in  the  resulting 
estimate.  This,  in  turn,  would  ensure  the  positive  definiteness  of  the  matrix 
in  question.  This  reasoning  can  be  extended  inductively  for  all  i  by  observing 
from  (67)  that  &n/n  j  must  be  positive  definite  if  Pn  j  has  that  property  (it 
is  implicitly  assumed  that  [9gn_j/3xn  j]  is  invertible).  Thus,  in  practical 
situations  where  weighting  matrices  have  been  sensibly  chosen,  a  matrix 
singularity  should  not  arise. 

An  interesting  interpretation  of  P^  ^  and  P^  can  be  obtained  by  com¬ 
paring  (64)  and  (65)  with  (60)  and  (62).  These  relations  can  readily  be  used  to 
show  that  P.y.  j  is  the  appropriate  "a  priori"  weighting  matrix  which  trans¬ 
forms  the  basic  problem  of  (30)  into  a  sequence  of  single  stage  problems. 

That  is,  minimizing 


vilyi  -*i)  pi/i-ilyi •  xi]  +  rlzi  •WlX'i'i'W1  (79) 


for  i  =  2,3,  .  .  .  ,  n  is  equivalent  to  minimizing  (30).  In  this  context,  it  can 
also  be  shown  that  P.  =  3(3J.  /  3x. )/ 3x. . 

l  ill 

The  formula  (79)  leads  to  a  heuristic  discussion  of  the  behavior  of 
Pi/i  j .  In  well  behaved  filtering  problems  where  large  amounts  of  data  are 
available  it  seems  that  there  should  be  a  point  beyond  which  new  data  is 
redundant.  By  (79)  this  could  only  occur  if  P^  j  reached  a  steady-state 
value  which  is  small"  compared  to  R^.  Certainly  if  I\  jQ^  j I\  j  is  large, 
this  could  not  happen.  Note,  however,  that  new  data  should  never  be  redundant 
in  such  cases  since  Q  is  chosen  to  be  large  when  large  disturbances  are  antici¬ 
pated  or  when  there  is  little  confidence  in  the  system  dynamic  model.  In 
general,  the  behavior  of  the  matrix  P.y^  i  *n  practical  problems  should  pro¬ 
vide  some  indication  as  to  the  adequacy  of  the  filtering  process.  Note  that 
establishing  conditions  for  convergence  of  the  sequential  estimation  formulas 
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of  Table  II  is  an  open  problem.  It  is  conjectured  that  the  solution  is  the  key 
to  nonlinear  observability. 

The  modeling  problem  will  not  be  discussed  in  detail  because  the  com¬ 
ments  made  for  continuous  systems  in  that  regard  apply  here  essentially 
unchanged.  Note  that  the  appropriate  filtered  and  smoothed  estimates  of 
the  dynamic  error  are,  respectively, 


a  Sh-  1 

r.  .w.  ,  =  P.  — —  Rr  1  [ z.  -  h.(x. )] 
l-l  l-l  1  dx  1  1  1  i  i 

l 


A  * 

r.  .w.  (  =  r.  .q.  .r.  a. 

i-1  i-1/n  l-l  l-l  l-l  l 


Comparison  With  Differential  Corrections 


(2  £  i  <  n)  (81) 


Adequate  results  have  often  been  obtained  in  practical  nonlinear  estima¬ 
tion  problems  using  linear  techniques  such  as  the  Kalman  Filter  in  conjunction 
with  the  method  of  differential  corrections  [7,8]  .  The  basic  philosophy  of  this 
approach  can  be  easily  demonstrated  and  compared  to  the  method  proposed 
here  using  a  single  stage  estimation  problem  as  an  example. 

To  minimize 

J  =  j  [z  -  h(y)J  R_1[z  -  h(y)]  (82) 

by  the  method  of  differential  corrections,  it  is  assumed  that  the  minimizing 
value  of  y  is  close  to  some  nominal  x  so  that  h(y)  can  be  adequately  approxi¬ 
mated  as  a  linear  function  of  y-  x.  That  is, 

J  »  jjz  -  h(x)  -  ~  [  y  -  x]  j  R"1  £z  -  h(x)  -  —  [  y  -  x]J  (83) 


Then,  3(J)/3y  =  0  gives  the  estimate 


[2  -  h(x)] 


(84) 


Note  that  the  linearizing  approximation  is  made  before  minimization  is 
attempted. 

The  alternative  presented  here  is  essentially  one  of  first  minimizing 
and  then  approximating.  Thus,  referring  to  (82), 


R‘l[z  -  h(y ) ] 


=  0 


(85) 


Now  assume  that  3(J)/3y  is  a  linear  function  of  y  -  x  and  obtain 


.  -i- 

(—■  R"  l[z  -  h(x)A 

~~  R_1[z  -  h(x)] 

(86) 

3x 

'  3x 

3x 

or 

a  —  .  3h '  „-l  3h  3  /  9h  \  _  - 1  r  L(— .  1 1  3h  _  - 1 .  . 

x  =  x  +  I  -3-  R  — —  (  -3-)  R  [z  -  h(x) ]  R  [z  -  h(x)]  (87) 

(_  3x  3x  3x  \  3x  /  J  3x 

Compare  (84)  and  (87)  and  note  the  term  missing  in  (84).  The  basic  result 
(87)  can  also  be  derived  by  assuming  J  to  be  a  quadratic  function  of  y  -  x 
and  then  minimizing  that  approximation. 

The  observations  made  above  apply  in  principle  to  the  much  more 
complex  dynamic  estimation  problem.  Thus,  use  of  the  formulas  presented 
in  this  paper  instead  of  those  for  the  differentially  corrected  Kalman  filter 
should  provide  imprr  2d  accuracy  and  computational  efficiency  in  practical 
problems . 
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Unusual  Cases 


The  direct  use  of  the  basic  nonlinear  filtering  algorithm  of  (63)-(67) 
should  provide  efficient  and  accurate  state  estimation  for  most  practical  prob¬ 
lems. 

Some  difficulty  could  arise,  however,  when  the  data  taken  in  the  initial 
stages  are  of  a  particularly  poor  quality.  The  following  heuristic  "start  up" 
procedure  is  recommended  in  such  cases. 

Apply  (63)-(67)  for  several  stages.  Use  the  smoothing  algorithm  (68)- 
(70)  to  obtain  Xj.  If  =  0,  continue  filtering  from  the  point  it  was  interrupted. 
On  the  other  hand,  if  X j  ^  0,  two  choices  remain:  either  continue  with  filter¬ 
ing  for  a  few  more  stages  and  smooth  again,  or  restart  the  problem  using  the 
smoothed  result  at  the  first  stage  for  linearization.  That  is,  let  [8] 


pi  = 


i  -  p 


1/0  £ 


l/n  '  1/n 


Pi/0  (88) 


A  A 
Xj  =  X 


l/n+Pl 


a  Ri  [*i  *  hi(xi/n>]  +  Pi/0^x0  ‘  Xl/J 


L^l/n 


(89) 


Continue  this  "looping"  procedure  until  X^  =  0. 

Another  special  case  of  interest  is  encountered  if  large  amounts  of 
reliable  data  are  available.  In  such  situations  the  filtering  algorithm  of 
Table  II  may  "converge"  sufficiently  at  some  point  so  that  the  term 
[  0(flhi/8Xi)/8Xi]R7*[z^  -  h^x^)]  can  be  omitted  and  the  formulas  (74),  (76),  and 
(77)  used  in  place  of  (64)  and  (65).  Of  course,  this  would  only  be  advantageous 
if  the  dimension  of  the  measurement  vector  q  is  less  them  that  of  the  state 
vector  p. 


-24- 


•  fl 


APPENDIX  A.  CONTINUOUS  LINEAR  ESTIMATION 

Although  the  solutions  to  linear  estimation  problems  are  well  known[  2] , 
they  are  derived  by  a  special  technique  here  and  in  Appendix  B  specifically 
to  provide  insight  and  motivation  for  the  techniques  used  to  solve  the  cor¬ 
responding  nonlinear  estimation  problems. 

The  linear  version  of  the  continuous  least  squares  estimation  problem 
is  obtained  if  in  (l)-(2)  f(x,t)  =  F(t)x  and  h(x,  t)  =  H(t)x.  Then  (5)-(7)  become 

y  =  Fy  -  rQr*\  (A-  1) 

X  =  -F*X  +  H*R_1[z  -  Hy]  (A-2) 

\(t0)  =  -Pq1  [y (tQ)  -  xQ]  (X(T)  =  0)  (A-3) 


By  expressing  (A-l)  and  (A-2)  in  the  augmented  vector-matrix  form 
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(A-4) 


and  partitioning  the  corresponding  transition  matrix  accordingly,  the  solution 
can  be  expressed  as 


y(t) 


=  *lt(t't0)y(t0)  +  ®l2(t't0)Mtl)  +  /t*12( t,o-)HVlz(«r)  do-  (A-5) 

t  n 


Mt)  =  *->1<t»Vy(t0>  +  ®22(t,t0)X(t1)  +  J *  *22(t,<r)H*R’lz(<r) 

t  A 


d<r 


(A-6) 
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This  result  can  be  transformed  to  a  much  more  useful  form  as  follows.  Solve 
(A-6)  for  y(tQ)  and  substitute  the  results  in  (A- 5).  Then 

y(t)  =  j(t,t0)*21(t,t0)\(t)  +  2^**0^  -  *11^' *0^*21  ^t,t0^®22^t,t0^^t0^ 

t  t 

-®11(t,t0)*211(t,t0)/’  ®22(t,(r)H*R‘1z((r)d<r  +f  •12(ti«r)H*R’1*(x)ds- 

**0  \ 

(A-7) 

Since  (A-7)  holds  for  arbitrary  tg,  the  solution  to  equations  (A-l)  and  (A-2) 
for  the  boundary  condition  X(T)  =  0  can  be  expressed  in  the  general  form 


y(t)  =  -P(t)X(t)  +  lit) 


(A-8) 


with  P(t)  and  g(t)  as  yet  unspecified.  Differentiating  (A-8)  and  substituting 
from  (A-l)  and  (A-2)  yields 

-g  +  FI  -  PH*R_1Hg  +  PH*R_1z  =  [-P  +  FP  +  PF*  -  PH*R_1HP  +  rQr*]X 

(A-9) 


For  this  to  hold  for  all  X  requires 


g  =  Fg  +  PH*R-1[z  -  Hg] 


(A- 10) 


P  =  FP  +  PF*  +  rQT*  -  PH*R-1HP 


(A-l  1) 


From  (A-8)  and  (A-3)  the  initial  conditions  for  these  equations  are 


seen  to  be 


.  ,  .  A 

g(tn)  =  Xn 


(A-12) 


P‘V  ■  P0 


(A- 13) 


Further,  since  X(T)  =  0,  this  gives 

y(T)  =  4<T)  (A- 14) 

The  proper  interpretation  of  the  result  is  the  key  to  the  solution  desired. 
Recall  that  y  is  the  smoothed  estimate  in  the  interval  t^i  t  £  T.  From  this  it 
should  be  clear  that  £(t)  for  arbitrary  t  is  an  estimate  of  the  state  vector  based 
on  the  data  gathered  from  tg  to  t.  By  definition,  then,  this  is  the  "filtered" 
estimate  of  the  state.  Note  also  that  for  any  T,  the  term  £(T)  is  the  missing 
boundary  condition  on  (A-l)  so  that  a  backward  solution  of  (A-l)-(A-2)  using 

X(T)  =  0  and  y(T)  =  $(T)  (A-l  5) 

immediately  gives  the  complete  smoothing  solution  in  the  interval  t^  £  t  £T. 

Note  that  for  the  standard  statistical  model  where  in  (1)  and  (2)  w  ->nd  v 
are  white  noise  processes  with  spectral  density  matrices  Q  and  R,  respectively, 
the  filter  (A- 1 0)-(A- 1 1)  is  the  Kalman  filter.  In  any  case  it  is  the  optimal 
least  squares  filter. 

A  useful  interpretation  of  the  matrix  P  can  be  obtained  from  the  defining 
relationship  (A- 8).  Differentiating  with  respect  to  X  yields 


~(x(t/T)  -  x(t>)  =  -P(t) 


(A- 1 6) 


From  this  it  is  clear  that  a  "small"  P  indicates  that  the  difference  between 
the  smoothed  and  filtered  estimates  is  relatively  unaffected  by  changes  in  X. 
This  would  be  expected  to  occur  after  large  amounts  of  data  have  been  pro¬ 
cessed  or  when  the  a  priori  guess  of  the  state  vector  is  especially  consistent 
with  the  observed  data.  (Observability  is  implicitly  assumed  here. )  In  such 
cases,  the  smoothing  algorithm  (A-l),  (A-2)  and  (A- 15)  would  be  of  little  value. 
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On  the  other  hand,  if  P  were  "large",  (A-16)  clearly  indicates  an  improvement 
in  estimating  the  state  from  smoothing.  In  a  qualitative  sense,  then,  P  is  a 
measure  of  the  effectiveness  of  filtering  where  the  ultimate  standard  of 
"goodness"  is  the  smoothed  estimate. 

An  interesting  observation  can  be  made  from  (A-8).  The  main  motiva¬ 
tion  for  its  introduction  was  to  solve  the  filtering  problem.  If,  instead,  the 
object  had  been  to  solve  the  smoothing  (tracking)  problem  directly,  the  desired 
substitution  would  be 


Mt)  =  -S(t)y(t)  +  Ut)  (A- 17) 

This  gives  "feedback  control"  directly  when  substituted  in  (A-l).  Now,  com¬ 
paring  (A- 17)  with  (A-8)  it  is  immediately  obvious  that 

S  =  P"1  (A- 18) 

and,  indeed,  the  resulting  differential  equation  for  S  is  precisely  the  equation 
for  P  * .  This  is,  of  course,  the  well-known  duality  between  the  optimal 
control  and  optimal  filtering  problems.  Note  that  substitutions  of  the  form 
(A-17)  have  been  used  effectively  in  solving  nonlinear  TPBY  problems  by 
McReynolds  [  10  ] . 


APPENDIX  B.  DISCRETE  LINEAR  ESTIMATION 


The  linear  version  of  the  discrete  least  squares  estimation  problem  is 
obtained  if  in  (2 8 ) - ( 29)  g^(x.)  =  and  h^(x^)  =  ILx^  In  that  case  (34)-(37) 
becomes 


Xl  =®lX2  +  HlRl‘lzl  -  Hl',l1  'Po'lyl  -*0]  <B-1) 


-28- 


=  «A.X.  +  h!  r7  *f  z.  -  Ey.j 
1  i  l+l  i  i  1  i  i'iJ 


(B-2) 


y.  =  4.  (y.  ,  +  r.  .Q.  .  r.  A. 

ri  l-l  l-l  l-l  l-l  l-i  i 


(2  £  i  £  n) 


(B-3) 


X1  =  Xn+1  =  0 


(B-4) 


Solving  { B - 1 )  for  with  \j  =  0  immediately  gives 


Now  let 


v,  =  [i+pohi'r;,h,]''poKx2+po,!so+hir;,zi]  •  <b-5» 


pt  -[I  +  poHX'H,] 


(B-6i 


et  =[l  *  poH^‘Hi]'lpo[po'2o  +  HX‘2.] 

=  x0  +  [J  +  P0H1RI  Hi]  poh,»i‘[«|-hi: 


so  that  (B-5)  can  be  written  as 


(B-7) 


V1  =  p  1*1^2  +  ^1 


(B-81 


This  suggests  that  in  general  the  change  of  variable 


y.  =  P.$A.. ,  + 

i  ii  i+l 


( 1  £  i  £  n) 


(B-9) 
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H 


(where  and  are  arbitrary)  be  applied  to  (B-i)  and  (B-2).  Note  from 
(B-4),  however,  that  if  (B-9)  is  to  be  valid,  must  be  the  smoothed  state 
vector  at  the  end  of  an  i-stage  process.  By  definition  this  is  the  "filtered” 
state  vector  and,  therefore, 


4, 


(B-10) 


where  the  notation  x^  is  introduced  to  specifically  indicate  a  filtered  estimate. 
Note  that  for  i  =  n  =  1,  has  already  been  defined  correctly. 

To  apply  (B-9)  to  (B-2)  and  (B-3)  write 


_  *  .  .  A 

v.  .  =  P.  ,  ®.  .X.  +  x.  , 
*1-1  l-l  l-l  i  l-l 


(B- 1 1 ) 


*  a 

P.®.  X.,.  +  x. 
lii+i  l 


(B-12) 


Combining  this  with  (B-3)  gives 


ri.iQi.ir;>i  ■  *a\, + k  -  -  *i.A.t  (b-i3> 


By  first  applying  (B-2)  and  then  (B-12),  X^  and  xi  can  be  eliminated  from 
(B- 13)  to  give 


#  $  #  - 1 r  r  ❖  Ail 

r.  .Q.  ,r.  ,x...  +  r.  .q.  ,r.  .h. r.  z.  -  h.  p.®.x.. .  +  x. 
i-l  i-l  i-l  i+l  i-l  i-l  i-l  i  i  L  i  il  i  i  i+I  iJJ 


_  *  ,  A  .  A 

=  P.®.  X. , ,  +  x.  -  ®.  .x.  . 
l  l  i+l  l  i-l  i-l 


#  *  -1 

®.X...  +  H.  R.  1 

ii+l  li 


X  [zi  -  Hi  [pi*i  xi+ 1  +  *i]] 


(B- 14) 
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For  convenience  define 


A  *  * 

p.,.  .  =  r.  ,q.  .r.  .  +  *.  .p.  .4.  , 
x/i-l  1-1  1-1  1-1  1-1  1-1  1-1 


(B-15) 


Then  (B-14)  can  be  simplified  to 


❖  ❖  _  I  A  a 

P. ,.  +  P./.  .  H.  R.  z.  -  x.  +  4.  .x.  . 

l/i-l  l  l+l  l/i-l  ill  l  l-l  i-l 


sP.*fL,+P.,  .  H*R.~  *H.P.4f\..xl  +P.  ».  .H*R:1H.x.  (B-16) 

l  i  l+l  l/i-l  ii  ill  l+l  l/i-l  ii  ii 


For  i  =  n,  the  last  stage  of  an  estimation  problem,  X.^+^  =  0  by  (B-4)  so  that 
(B-16)  becomes 


A 

X 


n 


_A 

4  .x  .  + 
n-1  n-1 


I  +  P  ,  .  H*R-1H  1  1P  ,  .  H  '  R" 1  [  z  -  H  4  ,  x  .] 

n/n-1  n  n  nj  n/n-1  n  n  1  n  n  n-1  n-1 J 


(B- 17) 


Suppose  now  that  additional  data  is  acquired,  i.e.,  the  last  stage 
becomes  n  +  1 .  Then  Xn+£  =  ®  ^n+i  ^  0*  so  that  equating  coefficients  of 
in  (B-16)  for  i  =  n  and  simplifying  gives 


P 

n 


❖  -  I 

I+P  .  .HR1 
n/n-1  n  n 


(B- 18) 


Thus  the  formulas  (B-17),  (B-18),  and  (B-15)  with  the  initial  conditions  (B-6) 
and  (B-7)  can  be  used  recursively  for  least  squares  filtering  of  an  arbitrary 
linear  process. 

Note  also  that  for  any  n,  the  term  xr  is  the  missing  boundary  condition 
on  (B-4)  so  that  a  backward  solution  of  (B-2)-(B-3)  using 

yn  =  £n  and  Xn+I=°  <B-l9) 

immediately  gives  the  complete  smoothing  solution. 
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APPENDIX  C .  CONTINUOUS  MEASUREMENT  LIMIT  FOR 
THE  DISCRETE  NONLINEAR  FILTER 


It  will  be  shown  here  that  as  the  time  between  successive  measurements 
approaches  zero,  the  discrete  filtering  formulas  of  Table  II  approach  the  con¬ 
tinuous  ones  of  Table  I. 

To  begin,  some  basic  relationships  must  be  established  between  the 
various  quantities  involved.  Integrating  (4)  between  and  t  with 
At  =  t.  -  t.  j  assumed  to  be  small  yields 

yi  yi-l  +  ^i-l^i-l^*  +  (C- 

Comparing  this  with  (31)  immediately  gives 


=  yi-l + 

(C-2) 

ri-lui-l  =  r(ti-l)ulti-l)At 

(C-3) 

whereupon 

^8i~  1  -  t  i  «  -  , 

8Vl  Syi-1 

(C-4) 

To  force  the  two  least  squares  criteria  (3)  and  (30) 

is  clear  that 

to  be  equivalent,  it 

R7l  =  R_1(t.)At 

1  1 

(C-5) 

-i  Q'1(V 

i- 1  At 

(C-6) 
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where  (C-3)  was  used  to  establish  the  last  result.  Combining  (C-4),  (C-6), 
and  (67)  gives 


pi/i-i  =  r(ti-i)Q<ti-i>r '<ti.i)At+ pi-i  +  a^7Pi-iAt+pi.i  W 


9f 


8f 


At 


i-1 


+  -“L.  P. 


8*i-l  i_1  3yi-l 


8f  A*2 

At 


(C-7) 


In  order  to  use  this  relationship  effectively,  (65)  must  be  expressed  in  the 
form 


=  P 


i/i-1 


(C-8) 


where  the  notational  abbreviation  (46)  has  been  used  for  convenience.  Now 
use  (C-7)  to  eliminate  the  first  term  on  the  right  of  (C-8).  Then 


P.  -  P.  . 

i  l-l 


df 


df 


St-1-  ■  8^7  Pi-1  +  Pi-1  8^  +  «H.|» 


df 


8Vi-l  U1  ®Vi-l 


At  +  P. 


‘8r|.  ,  1  8r|.  1-1 

1  ±_  T  p  1  p 

i/i_i  8x  At  ’  i/i’i  £  {/i‘ 

L  i  J  *•  i-J 


(C-9) 
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